HS project

IBL and ORS signature

This file is used to identify specific markers for IBL and ORS.

library(dplyr)
library(patchwork)
library(ggplot2)
library(ComplexHeatmap)

.libPaths()

## [1] "/usr/local/lib/R/library"

Preparation

In this section, we set the global settings of the analysis. We will store data there :

out_dir = "."

We load the dataset :

sobj = readRDS(paste0(out_dir, "/hs_hd_sobj.rds"))
sobj

## An object of class Seurat 
## 20003 features across 12111 samples within 1 assay 
## Active assay: RNA (20003 features, 2000 variable features)
##  6 dimensional reductions calculated: RNA_pca, RNA_pca_38_tsne, RNA_pca_38_umap, harmony, harmony_38_umap, harmony_38_tsne

We load the sample information :

sample_info = readRDS(paste0(out_dir, "/../1_metadata/hs_hd_sample_info.rds"))
project_names_oi = sample_info$project_name

graphics::pie(rep(1, nrow(sample_info)),
              col = sample_info$color,
              labels = sample_info$project_name)

Here are custom colors for each cell type :

color_markers = readRDS(paste0(out_dir, "/../1_metadata/hs_hd_color_markers.rds"))

data.frame(cell_type = names(color_markers),
           color = unlist(color_markers)) %>%
  ggplot2::ggplot(., aes(x = cell_type, y = 0, fill = cell_type)) +
  ggplot2::geom_point(pch = 21, size = 5) +
  ggplot2::scale_fill_manual(values = unlist(color_markers), breaks = names(color_markers)) +
  ggplot2::theme_classic() +
  ggplot2::theme(legend.position = "none",
                 axis.line = element_blank(),
                 axis.title = element_blank(),
                 axis.ticks = element_blank(),
                 axis.text.y = element_blank(),
                 axis.text.x = element_text(angle = 30, hjust = 1))

This is the projection of interest :

name2D = "harmony_38_tsne"

We design a custom function to make a histogram and a wordcloud to visualize differentially expressed genes :

hist_wc_fun = function(mark, col) {
  cut_colors = c("firebrick4", "firebrick2", "indianred1", "darksalmon",
                 "lightpink", "gray50", "khaki3", "darkolivegreen1",
                 "olivedrab1", "chartreuse2", "chartreuse4")
  cut_colors_cont = c(rev(RColorBrewer::brewer.pal(name = "Reds", n = 9)[c(5:9)]),
                      RColorBrewer::brewer.pal(name = "Greens", n = 9)[c(5:9)])
  
  if (col == "pct.1") {
    mark$pct = mark$pct.1
    mark$pct_cut = mark$pct.1_cut
  } else if (col == "pct.2") {
    mark$pct = mark$pct.2
    mark$pct_cut = mark$pct.2_cut
    cut_colors = rev(cut_colors)
    cut_colors_cont = rev(cut_colors_cont)
  } else {
    stop("col must be either pct.1 or pct.2")
  }
  
  p_hist = ggplot2::ggplot(mark, mapping = aes(x = avg_logFC, fill = pct_cut)) +
    ggplot2::geom_histogram(binwidth = 0.05) +
    ggplot2::scale_fill_manual(breaks = levels(mark$pct_cut),
                               values = cut_colors,
                               name = paste0(col, "_cut")) +
    ggplot2::theme_classic()
  
  p_wc = ggplot2::ggplot(mark, aes(label = gene_name, size = avg_logFC, color = pct)) +
    ggwordcloud::geom_text_wordcloud_area(show.legend = TRUE, seed = 1) +
    ggplot2::scale_color_gradientn(colors = cut_colors_cont,
                                   name = col) +
    ggplot2::scale_size_area(max_size = 5) +
    ggplot2::theme_minimal() +
    ggplot2::guides(size = "none")
  
  p = patchwork::wrap_plots(p_hist, p_wc, nrow = 1)
  
  return(p)
}

Visualization

Gene expression

We visualize gene expression for some markers :

features = c("percent.mt", "percent.rb", "nFeature_RNA")

plot_list = lapply(features, FUN = function(one_gene) {
  Seurat::FeaturePlot(sobj, features = one_gene,
                      reduction = name2D) +
    ggplot2::theme(aspect.ratio = 1) +
    ggplot2::scale_color_gradientn(colors = aquarius::color_gene) +
    Seurat::NoAxes()
})

patchwork::wrap_plots(plot_list, ncol = 3)

Cluster type Clusters and cell type

We visualize clusters and cell type :

cluster_plot = Seurat::DimPlot(sobj, group.by = "seurat_clusters",
                               reduction = name2D, label = TRUE) +
  ggplot2::labs(title = "Cluster ID") +
  Seurat::NoAxes() +
  ggplot2::theme(aspect.ratio = 1,
                 plot.title = element_text(hjust = 0.5))

cell_type_plot = Seurat::DimPlot(sobj, group.by = "cell_type",
                                 reduction = name2D, label = FALSE) +
  ggplot2::scale_color_manual(values = color_markers,
                              breaks = names(color_markers)) +
  ggplot2::labs(title = "Cell type") +
  Seurat::NoAxes() +
  ggplot2::theme(aspect.ratio = 1,
                 plot.title = element_text(hjust = 0.5))

cell_type_plot | cluster_plot

We summarize major cell type by cluster :

cell_type_clusters = sobj@meta.data[, c("cell_type", "seurat_clusters")] %>%
  table() %>%
  prop.table(., margin = 2) %>%
  apply(., 2, which.max)
cell_type_clusters = setNames(levels(sobj$cell_type)[cell_type_clusters],
                              nm = names(cell_type_clusters))

We define cluster type :

sobj$cluster_type = cell_type_clusters[sobj$seurat_clusters] %>%
  as.factor()
table(sobj$cluster_type, sobj$cell_type)

##                   
##                    CD4 T cells CD8 T cells Langerhans cells macrophages B cells
##   B cells                    0           0                0           0      37
##   CD4 T cells              774          52                2           1       3
##   CD8 T cells               86         545                1           0       4
##   HFSC                       0           0                0           0       2
##   IBL                        0           0                0           0       1
##   IFE                        0           0                0           0       1
##   IRS                        0           0                0           0       0
##   Langerhans cells          15           0              253          51       1
##   ORS                        2           0                0           0       3
##   cortex                     0           0                0           0       0
##   cuticle                    1           0                0           0       0
##   macrophages                9           0               26         423       1
##   medulla                    0           0                0           0       1
##   proliferative              0           0                3           0      12
##   sebocytes                  0           0                0           0       0
##                   
##                    cuticle cortex medulla  IRS proliferative  IBL  ORS  IFE
##   B cells                0      0       1    0             0    0    0    0
##   CD4 T cells            2      0       3    3             2    1    0    0
##   CD8 T cells            0      0       0    0             0    0    0    0
##   HFSC                   4      0       2    2             9   41  104    3
##   IBL                    4      0       5    1            23 1502   11  113
##   IFE                    1      0       0    0             2   25   66  498
##   IRS                    0      0       0  155             4    0    0    0
##   Langerhans cells       1      1       0    2            23    2    0    0
##   ORS                    5      2      15    6            23   21 1779   57
##   cortex               105    199       4    0             0    0    0    0
##   cuticle              804     29      21    0             1    0    0    0
##   macrophages            0      0       0    1             0    0    0    0
##   medulla               26      1     335    0            17    0    0    0
##   proliferative        318      8      42  126          1384   44   75   71
##   sebocytes              0      1       1    0             1    0   63    1
##                   
##                    HFSC sebocytes
##   B cells             0         0
##   CD4 T cells         0         0
##   CD8 T cells         0         0
##   HFSC             1286         1
##   IBL                23         7
##   IFE                 7        13
##   IRS                 0         0
##   Langerhans cells    1         2
##   ORS                44        16
##   cortex              0         0
##   cuticle             0         2
##   macrophages         0         0
##   medulla             1         0
##   proliferative      24        18
##   sebocytes           0       154

We compare cluster annotation and cell type annotation :

cell_type_plot

p2 = Seurat::DimPlot(sobj, group.by = "cluster_type",
                     reduction = name2D, cols = color_markers) +
  ggplot2::labs(title = "Cluster type") +
  Seurat::NoAxes() +
  ggplot2::theme(aspect.ratio = 1,
                 plot.title = element_text(hjust = 0.5))

patchwork::wrap_plots(cell_type_plot, p2, guides = "collect")

There is mis-annotation in ORS (small part close to IBL) and in IBL (small part in IFE), so we keep the single-cell level cell type annotation.

Differential expression

In this section, we perform DE between inner bulge layer (IBL) or outer root sheath (ORS), and all remaining cells. We save the results in a list :

list_results = list()

We change cell identities to cell type :

Seurat::Idents(sobj) = sobj$cell_type

table(Seurat::Idents(sobj))

## 
##      CD4 T cells      CD8 T cells Langerhans cells      macrophages 
##              887              597              285              475 
##          B cells          cuticle           cortex          medulla 
##               66             1270              241              429 
##              IRS    proliferative              IBL              ORS 
##              296             1489             1636             2098 
##              IFE             HFSC        sebocytes 
##              743             1386              213

IBL markers

In this section, we compared IBL to other cells.

group_name = "IBL_vs_all"

aquarius::plot_red_and_blue(sobj,
                            group1 = "IBL",
                            reduction = name2D) +
  ggplot2::labs(title = group_name)

Differential expression

We identify specific markers for each population :

mark = Seurat::FindMarkers(sobj, ident.1 = "IBL")

mark = mark %>%
  dplyr::filter(p_val_adj < 0.05) %>%
  dplyr::arrange(-avg_logFC, pct.1 - pct.2)
mark$gene_name = rownames(mark)

list_results[[group_name]]$mark = mark

dim(mark)

## [1] 970   6

head(mark, n = 20)

##                  p_val avg_logFC pct.1 pct.2     p_val_adj gene_name
## KRT16     0.000000e+00  3.892906 0.905 0.155  0.000000e+00     KRT16
## KRT6B     0.000000e+00  3.166639 0.921 0.275  0.000000e+00     KRT6B
## KRT17     0.000000e+00  2.954873 0.980 0.529  0.000000e+00     KRT17
## KRT6A     0.000000e+00  2.860010 0.752 0.260  0.000000e+00     KRT6A
## S100A2    0.000000e+00  2.635463 0.992 0.789  0.000000e+00    S100A2
## KRT6C     0.000000e+00  2.494914 0.548 0.068  0.000000e+00     KRT6C
## KRT5      0.000000e+00  1.988478 0.993 0.749  0.000000e+00      KRT5
## FABP5     0.000000e+00  1.975885 0.989 0.668  0.000000e+00     FABP5
## GJB6      0.000000e+00  1.912305 0.908 0.525  0.000000e+00      GJB6
## KRT14     0.000000e+00  1.758304 0.992 0.699  0.000000e+00     KRT14
## NDUFA4L2  0.000000e+00  1.721714 0.680 0.212  0.000000e+00  NDUFA4L2
## CST6      0.000000e+00  1.666433 0.384 0.076  0.000000e+00      CST6
## GJB2      0.000000e+00  1.640163 0.826 0.503  0.000000e+00      GJB2
## HSPB1     0.000000e+00  1.630779 0.985 0.835  0.000000e+00     HSPB1
## TM4SF1    0.000000e+00  1.597284 0.801 0.243  0.000000e+00    TM4SF1
## CALML3    0.000000e+00  1.590650 0.958 0.685  0.000000e+00    CALML3
## SDC1      0.000000e+00  1.582781 0.829 0.544  0.000000e+00      SDC1
## S100A16   0.000000e+00  1.547152 0.955 0.720  0.000000e+00   S100A16
## ANXA2     0.000000e+00  1.518822 0.991 0.852  0.000000e+00     ANXA2
## LGALS7B  3.086968e-238  1.508193 0.748 0.465 6.174861e-234   LGALS7B

How many genes enriched in IBL ?

mark_IBL = mark %>%
  dplyr::filter(avg_logFC > 0)

nrow(mark_IBL)

## [1] 336

We cut pct.1 and pct.2 by bins of 0.1 :

mark_IBL$pct.1_cut = cut(mark_IBL$pct.1, breaks = 10)
mark_IBL$pct.2_cut = cut(mark_IBL$pct.2, breaks = 10)

head(mark_IBL)

##        p_val avg_logFC pct.1 pct.2 p_val_adj gene_name     pct.1_cut
## KRT16      0  3.892906 0.905 0.155         0     KRT16 (0.905,0.995]
## KRT6B      0  3.166639 0.921 0.275         0     KRT6B (0.905,0.995]
## KRT17      0  2.954873 0.980 0.529         0     KRT17 (0.905,0.995]
## KRT6A      0  2.860010 0.752 0.260         0     KRT6A (0.726,0.815]
## S100A2     0  2.635463 0.992 0.789         0    S100A2 (0.905,0.995]
## KRT6C      0  2.494914 0.548 0.068         0     KRT6C (0.547,0.637]
##              pct.2_cut
## KRT16      (0.102,0.2]
## KRT6B      (0.2,0.297]
## KRT17    (0.493,0.591]
## KRT6A      (0.2,0.297]
## S100A2   (0.786,0.884]
## KRT6C  (0.00302,0.102]

Visualization

We make a histogram for pct.1, pct.2 and avg_logFC.

hist_wc_fun(mark_IBL, "pct.1")

hist_wc_fun(mark_IBL, "pct.2")

Selection

The best markers have high pct.1 and low pct.2 :

mark_IBL = mark_IBL %>%
  dplyr::filter(pct.1 > 0.6 & pct.2 < 0.3)

list_results[[group_name]]$choosen_ones = mark_IBL

mark_IBL

##          p_val avg_logFC pct.1 pct.2 p_val_adj gene_name     pct.1_cut
## KRT16        0  3.892906 0.905 0.155         0     KRT16 (0.905,0.995]
## KRT6B        0  3.166639 0.921 0.275         0     KRT6B (0.905,0.995]
## KRT6A        0  2.860010 0.752 0.260         0     KRT6A (0.726,0.815]
## NDUFA4L2     0  1.721714 0.680 0.212         0  NDUFA4L2 (0.637,0.726]
## TM4SF1       0  1.597284 0.801 0.243         0    TM4SF1 (0.726,0.815]
## LYPD3        0  1.408274 0.732 0.289         0     LYPD3 (0.726,0.815]
##            pct.2_cut
## KRT16    (0.102,0.2]
## KRT6B    (0.2,0.297]
## KRT6A    (0.2,0.297]
## NDUFA4L2 (0.2,0.297]
## TM4SF1   (0.2,0.297]
## LYPD3    (0.2,0.297]

We visualize expression levels of those genes on the projection :

plot_list = lapply(rownames(mark_IBL), FUN = function(one_gene) {
  Seurat::FeaturePlot(sobj, features = one_gene,
                      reduction = name2D) +
    ggplot2::theme(aspect.ratio = 1) +
    ggplot2::scale_color_gradientn(colors = aquarius::color_gene) +
    Seurat::NoAxes()
})

patchwork::wrap_plots(plot_list, ncol = 3)

ORS markers

In this section, we compared ORS to other cells.

group_name = "ORS_vs_all"

aquarius::plot_red_and_blue(sobj,
                            group1 = "ORS",
                            reduction = name2D) +
  ggplot2::labs(title = group_name)

Differential expression

We identify specific markers for each population :

mark = Seurat::FindMarkers(sobj, ident.1 = "ORS")

mark = mark %>%
  dplyr::filter(p_val_adj < 0.05) %>%
  dplyr::arrange(-avg_logFC, pct.1 - pct.2)
mark$gene_name = rownames(mark)

list_results[[group_name]]$mark = mark

dim(mark)

## [1] 739   6

head(mark, n = 20)

##                 p_val avg_logFC pct.1 pct.2     p_val_adj gene_name
## KRT15    0.000000e+00 1.2525184 0.929 0.367  0.000000e+00     KRT15
## ATP1B3   0.000000e+00 1.2365428 0.973 0.628  0.000000e+00    ATP1B3
## CCL2     0.000000e+00 1.1973508 0.506 0.056  0.000000e+00      CCL2
## DST      0.000000e+00 1.1908763 0.981 0.352  0.000000e+00       DST
## IMPA2    0.000000e+00 1.0522498 0.924 0.455  0.000000e+00     IMPA2
## COL17A1  0.000000e+00 0.9568857 0.923 0.299  0.000000e+00   COL17A1
## TXNIP    0.000000e+00 0.9062837 0.955 0.547  0.000000e+00     TXNIP
## S100A9   0.000000e+00 0.9041040 0.691 0.231  0.000000e+00    S100A9
## ZFP36L2  0.000000e+00 0.8729918 0.968 0.655  0.000000e+00   ZFP36L2
## IFITM3   0.000000e+00 0.8616994 0.909 0.433  0.000000e+00    IFITM3
## NEAT1    0.000000e+00 0.8281953 0.956 0.672  0.000000e+00     NEAT1
## MOXD1    0.000000e+00 0.8157750 0.764 0.126  0.000000e+00     MOXD1
## CEBPD    0.000000e+00 0.8026200 0.886 0.439  0.000000e+00     CEBPD
## AQP3     0.000000e+00 0.7900901 0.796 0.263  0.000000e+00      AQP3
## PTN     5.159311e-271 0.7884902 0.567 0.234 1.032017e-266       PTN
## ALDH3A1  0.000000e+00 0.7620031 0.607 0.087  0.000000e+00   ALDH3A1
## ARL4A    0.000000e+00 0.7563103 0.874 0.504  0.000000e+00     ARL4A
## GPX2     0.000000e+00 0.7468617 0.633 0.112  0.000000e+00      GPX2
## LIMA1    0.000000e+00 0.7441636 0.893 0.379  0.000000e+00     LIMA1
## FST      0.000000e+00 0.7440017 0.513 0.064  0.000000e+00       FST

How many genes enriched in ORS ?

mark_ORS = mark %>%
  dplyr::filter(avg_logFC > 0)

nrow(mark_ORS)

## [1] 275

We cut pct.1 and pct.2 by bins of 0.1 :

mark_ORS$pct.1_cut = cut(mark_ORS$pct.1, breaks = 10)
mark_ORS$pct.2_cut = cut(mark_ORS$pct.2, breaks = 10)

head(mark_ORS)

##         p_val avg_logFC pct.1 pct.2 p_val_adj gene_name     pct.1_cut
## KRT15       0 1.2525184 0.929 0.367         0     KRT15 (0.911,0.999]
## ATP1B3      0 1.2365428 0.973 0.628         0    ATP1B3 (0.911,0.999]
## CCL2        0 1.1973508 0.506 0.056         0      CCL2 (0.479,0.566]
## DST         0 1.1908763 0.981 0.352         0       DST (0.911,0.999]
## IMPA2       0 1.0522498 0.924 0.455         0     IMPA2 (0.911,0.999]
## COL17A1     0 0.9568857 0.923 0.299         0   COL17A1 (0.911,0.999]
##             pct.2_cut
## KRT15   (0.309,0.407]
## ATP1B3  (0.603,0.701]
## CCL2    (0.013,0.112]
## DST     (0.309,0.407]
## IMPA2   (0.407,0.505]
## COL17A1  (0.21,0.309]

Visualization

We make a histogram for pct.1, pct.2 and avg_logFC.

hist_wc_fun(mark_ORS, "pct.1")

hist_wc_fun(mark_ORS, "pct.2")

Selection

The best markers have high pct.1 and low pct.2 :

mark_ORS = mark_ORS %>%
  dplyr::filter(pct.1 > 0.6 & pct.2 < 0.2)

list_results[[group_name]]$choosen_ones = mark_ORS

mark_ORS

##         p_val avg_logFC pct.1 pct.2 p_val_adj gene_name     pct.1_cut
## MOXD1       0 0.8157750 0.764 0.126         0     MOXD1 (0.738,0.825]
## ALDH3A1     0 0.7620031 0.607 0.087         0   ALDH3A1 (0.566,0.652]
## GPX2        0 0.7468617 0.633 0.112         0      GPX2 (0.566,0.652]
## AHNAK2      0 0.7113339 0.749 0.184         0    AHNAK2 (0.738,0.825]
## CDH13       0 0.6071029 0.741 0.158         0     CDH13 (0.738,0.825]
## LAMB3       0 0.6026591 0.740 0.163         0     LAMB3 (0.738,0.825]
## S100A8      0 0.5772416 0.659 0.196         0    S100A8 (0.652,0.738]
## NBL1        0 0.4801919 0.627 0.159         0      NBL1 (0.566,0.652]
## CCDC3       0 0.4355107 0.601 0.182         0     CCDC3 (0.566,0.652]
## SPARC       0 0.4180143 0.607 0.180         0     SPARC (0.566,0.652]
##             pct.2_cut
## MOXD1    (0.112,0.21]
## ALDH3A1 (0.013,0.112]
## GPX2    (0.013,0.112]
## AHNAK2   (0.112,0.21]
## CDH13    (0.112,0.21]
## LAMB3    (0.112,0.21]
## S100A8   (0.112,0.21]
## NBL1     (0.112,0.21]
## CCDC3    (0.112,0.21]
## SPARC    (0.112,0.21]

We visualize expression levels of those genes on the projection :

plot_list = lapply(rownames(mark_ORS), FUN = function(one_gene) {
  Seurat::FeaturePlot(sobj, features = one_gene,
                      reduction = name2D) +
    ggplot2::theme(aspect.ratio = 1) +
    ggplot2::scale_color_gradientn(colors = aquarius::color_gene) +
    Seurat::NoAxes()
})

patchwork::wrap_plots(plot_list, ncol = 3)

Save

We save the list of results :

saveRDS(list_results, file = paste0(out_dir, "/ibl_ors_markers.rds"))

We also save as XLSX file :

list_results2 = list(IBL_vs_all = list_results$IBL_vs_all$mark,
                     IBL_vs_all_selection = list_results$IBL_vs_all$choosen_ones,
                     ORS_vs_all = list_results$ORS_vs_all$mark,
                     ORS_vs_all_selection = list_results$ORS_vs_all$choosen_ones)

openxlsx::write.xlsx(list_results2, file = paste0(out_dir, "/ibl_vs_ors_markers.xlsx"))

R Session

show

## R version 3.6.3 (2020-02-29)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 20.04.6 LTS
## 
## Matrix products: default
## BLAS:   /usr/local/lib/R/lib/libRblas.so
## LAPACK: /usr/local/lib/R/lib/libRlapack.so
## 
## locale:
## [1] C
## 
## attached base packages:
## [1] grid      stats     graphics  grDevices utils     datasets  methods  
## [8] base     
## 
## other attached packages:
## [1] ComplexHeatmap_2.14.0 ggplot2_3.3.5         patchwork_1.1.2      
## [4] dplyr_1.0.7          
## 
## loaded via a namespace (and not attached):
##   [1] softImpute_1.4              graphlayouts_0.7.0         
##   [3] pbapply_1.4-2               lattice_0.20-41            
##   [5] haven_2.3.1                 vctrs_0.3.8                
##   [7] usethis_2.0.1               dynwrap_1.2.1              
##   [9] blob_1.2.1                  survival_3.2-13            
##  [11] prodlim_2019.11.13          dynutils_1.0.5             
##  [13] later_1.3.0                 DBI_1.1.1                  
##  [15] R.utils_2.11.0              SingleCellExperiment_1.8.0 
##  [17] rappdirs_0.3.3              uwot_0.1.8                 
##  [19] dqrng_0.2.1                 jpeg_0.1-8.1               
##  [21] zlibbioc_1.32.0             pspline_1.0-18             
##  [23] pcaMethods_1.78.0           mvtnorm_1.1-1              
##  [25] htmlwidgets_1.5.4           GlobalOptions_0.1.2        
##  [27] future_1.22.1               UpSetR_1.4.0               
##  [29] laeken_0.5.2                leiden_0.3.3               
##  [31] clustree_0.4.3              parallel_3.6.3             
##  [33] scater_1.14.6               irlba_2.3.3                
##  [35] markdown_1.1                DEoptimR_1.0-9             
##  [37] tidygraph_1.1.2             Rcpp_1.0.9                 
##  [39] readr_2.0.2                 KernSmooth_2.23-17         
##  [41] carrier_0.1.0               promises_1.1.0             
##  [43] gdata_2.18.0                DelayedArray_0.12.3        
##  [45] limma_3.42.2                graph_1.64.0               
##  [47] RcppParallel_5.1.4          Hmisc_4.4-0                
##  [49] fs_1.5.2                    RSpectra_0.16-0            
##  [51] fastmatch_1.1-0             ranger_0.12.1              
##  [53] digest_0.6.25               png_0.1-7                  
##  [55] sctransform_0.2.1           cowplot_1.0.0              
##  [57] DOSE_3.12.0                 here_1.0.1                 
##  [59] TInGa_0.0.0.9000            ggraph_2.0.3               
##  [61] pkgconfig_2.0.3             GO.db_3.10.0               
##  [63] DelayedMatrixStats_1.8.0    gower_0.2.1                
##  [65] ggbeeswarm_0.6.0            iterators_1.0.12           
##  [67] DropletUtils_1.6.1          reticulate_1.26            
##  [69] clusterProfiler_3.14.3      SummarizedExperiment_1.16.1
##  [71] circlize_0.4.15             beeswarm_0.4.0             
##  [73] GetoptLong_1.0.5            xfun_0.35                  
##  [75] bslib_0.3.1                 zoo_1.8-10                 
##  [77] tidyselect_1.1.0            reshape2_1.4.4             
##  [79] purrr_0.3.4                 ica_1.0-2                  
##  [81] pcaPP_1.9-73                viridisLite_0.3.0          
##  [83] rtracklayer_1.46.0          rlang_1.0.2                
##  [85] hexbin_1.28.1               jquerylib_0.1.4            
##  [87] dyneval_0.9.9               glue_1.4.2                 
##  [89] RColorBrewer_1.1-2          matrixStats_0.56.0         
##  [91] stringr_1.4.0               lava_1.6.7                 
##  [93] europepmc_0.3               DESeq2_1.26.0              
##  [95] recipes_0.1.17              labeling_0.3               
##  [97] httpuv_1.5.2                class_7.3-17               
##  [99] BiocNeighbors_1.4.2         DO.db_2.9                  
## [101] annotate_1.64.0             jsonlite_1.7.2             
## [103] XVector_0.26.0              bit_4.0.4                  
## [105] mime_0.9                    aquarius_0.1.5             
## [107] Rsamtools_2.2.3             gridExtra_2.3              
## [109] gplots_3.0.3                stringi_1.4.6              
## [111] processx_3.5.2              gsl_2.1-6                  
## [113] bitops_1.0-6                cli_3.0.1                  
## [115] batchelor_1.2.4             RSQLite_2.2.0              
## [117] randomForest_4.6-14         tidyr_1.1.4                
## [119] data.table_1.14.2           rstudioapi_0.13            
## [121] org.Mm.eg.db_3.10.0         GenomicAlignments_1.22.1   
## [123] nlme_3.1-147                qvalue_2.18.0              
## [125] scran_1.14.6                locfit_1.5-9.4             
## [127] scDblFinder_1.1.8           listenv_0.8.0              
## [129] ggthemes_4.2.4              gridGraphics_0.5-0         
## [131] R.oo_1.24.0                 dbplyr_1.4.4               
## [133] BiocGenerics_0.32.0         TTR_0.24.2                 
## [135] readxl_1.3.1                lifecycle_1.0.1            
## [137] timeDate_3043.102           ggpattern_0.3.1            
## [139] munsell_0.5.0               cellranger_1.1.0           
## [141] R.methodsS3_1.8.1           proxyC_0.1.5               
## [143] visNetwork_2.0.9            caTools_1.18.0             
## [145] codetools_0.2-16            ggwordcloud_0.5.0          
## [147] Biobase_2.46.0              GenomeInfoDb_1.22.1        
## [149] vipor_0.4.5                 lmtest_0.9-38              
## [151] msigdbr_7.5.1               htmlTable_1.13.3           
## [153] triebeard_0.3.0             lsei_1.2-0                 
## [155] xtable_1.8-4                ROCR_1.0-7                 
## [157] BiocManager_1.30.10         scatterplot3d_0.3-41       
## [159] abind_1.4-5                 farver_2.0.3               
## [161] parallelly_1.28.1           RANN_2.6.1                 
## [163] askpass_1.1                 GenomicRanges_1.38.0       
## [165] RcppAnnoy_0.0.16            tibble_3.1.5               
## [167] ggdendro_0.1-20             cluster_2.1.0              
## [169] future.apply_1.5.0          Seurat_3.1.5               
## [171] dendextend_1.15.1           Matrix_1.3-2               
## [173] ellipsis_0.3.2              prettyunits_1.1.1          
## [175] lubridate_1.7.9             ggridges_0.5.2             
## [177] igraph_1.2.5                RcppEigen_0.3.3.7.0        
## [179] fgsea_1.12.0                remotes_2.4.2              
## [181] scBFA_1.0.0                 destiny_3.0.1              
## [183] VIM_6.1.1                   testthat_3.1.0             
## [185] htmltools_0.5.2             BiocFileCache_1.10.2       
## [187] yaml_2.2.1                  utf8_1.1.4                 
## [189] plotly_4.9.2.1              XML_3.99-0.3               
## [191] ModelMetrics_1.2.2.2        e1071_1.7-3                
## [193] foreign_0.8-76              withr_2.5.0                
## [195] fitdistrplus_1.0-14         BiocParallel_1.20.1        
## [197] xgboost_1.4.1.1             bit64_4.0.5                
## [199] foreach_1.5.0               robustbase_0.93-9          
## [201] Biostrings_2.54.0           GOSemSim_2.13.1            
## [203] rsvd_1.0.3                  memoise_2.0.0              
## [205] evaluate_0.18               forcats_0.5.0              
## [207] rio_0.5.16                  geneplotter_1.64.0         
## [209] tzdb_0.1.2                  caret_6.0-86               
## [211] ps_1.6.0                    DiagrammeR_1.0.6.1         
## [213] curl_4.3                    fdrtool_1.2.15             
## [215] fansi_0.4.1                 highr_0.8                  
## [217] urltools_1.7.3              xts_0.12.1                 
## [219] GSEABase_1.48.0             acepack_1.4.1              
## [221] edgeR_3.28.1                checkmate_2.0.0            
## [223] scds_1.2.0                  cachem_1.0.6               
## [225] npsurv_0.4-0                babelgene_22.3             
## [227] rjson_0.2.20                openxlsx_4.1.5             
## [229] ggrepel_0.9.1               clue_0.3-60                
## [231] rprojroot_2.0.2             stabledist_0.7-1           
## [233] tools_3.6.3                 sass_0.4.0                 
## [235] nichenetr_1.1.1             magrittr_2.0.1             
## [237] RCurl_1.98-1.2              proxy_0.4-24               
## [239] car_3.0-11                  ape_5.3                    
## [241] ggplotify_0.0.5             xml2_1.3.2                 
## [243] httr_1.4.2                  assertthat_0.2.1           
## [245] rmarkdown_2.18              boot_1.3-25                
## [247] globals_0.14.0              R6_2.4.1                   
## [249] Rhdf5lib_1.8.0              nnet_7.3-14                
## [251] RcppHNSW_0.2.0              progress_1.2.2             
## [253] genefilter_1.68.0           statmod_1.4.34             
## [255] gtools_3.8.2                shape_1.4.6                
## [257] HDF5Array_1.14.4            BiocSingular_1.2.2         
## [259] rhdf5_2.30.1                splines_3.6.3              
## [261] AUCell_1.8.0                carData_3.0-4              
## [263] colorspace_1.4-1            generics_0.1.0             
## [265] stats4_3.6.3                base64enc_0.1-3            
## [267] dynfeature_1.0.0            smoother_1.1               
## [269] gridtext_0.1.1              pillar_1.6.3               
## [271] tweenr_1.0.1                sp_1.4-1                   
## [273] ggplot.multistats_1.0.0     rvcheck_0.1.8              
## [275] GenomeInfoDbData_1.2.2      plyr_1.8.6                 
## [277] gtable_0.3.0                zip_2.2.0                  
## [279] knitr_1.41                  latticeExtra_0.6-29        
## [281] biomaRt_2.42.1              IRanges_2.20.2             
## [283] fastmap_1.1.0               ADGofTest_0.3              
## [285] copula_1.0-0                doParallel_1.0.15          
## [287] AnnotationDbi_1.48.0        vcd_1.4-8                  
## [289] babelwhale_1.0.1            openssl_1.4.1              
## [291] scales_1.1.1                backports_1.2.1            
## [293] S4Vectors_0.24.4            ipred_0.9-12               
## [295] enrichplot_1.6.1            hms_1.1.1                  
## [297] ggforce_0.3.1               Rtsne_0.15                 
## [299] shiny_1.7.1                 numDeriv_2016.8-1.1        
## [301] polyclip_1.10-0             lazyeval_0.2.2             
## [303] Formula_1.2-3               tsne_0.1-3                 
## [305] crayon_1.3.4                MASS_7.3-54                
## [307] pROC_1.16.2                 viridis_0.5.1              
## [309] dynparam_1.0.0              rpart_4.1-15               
## [311] zinbwave_1.8.0              compiler_3.6.3             
## [313] ggtext_0.1.0